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_. Abstract 

c/2 . We investigate using molecular dynamics the effect of the structure of nanoconfinenient for 



liquids with water-like anomalies and liquid-liquid phase transition (LLPT). We find that if the 
confinement is in an ordered matrix of nanoparticles (NPs) the anomalies are preserved, although 
the LLPT shifts to lower temperatures, higher pressures and higher densities with respect to bulk. 
On the contrary, if the NPs matrix is disordered, we find a drastically different phase diagram: 



^N ' the LLPT occurs only in a reduced interval of densities and the anomalies are washed out. To 

>■ 

}H ' understand this effect we calculate the changes in the system at the microscopic level. In all the 

different confinements considered here we observe a dramatic increase of density of liquid near the 

r^ I confining NPs. In the disordered case the confinement induces larger heterogeneity in the local 

density, responsible for the weakening of the LLPT and the disappearance of anomalies. 

^^ ■ PACS numbers: 64.70.Ja,65.20.-w, 66.10.C- 



Many experiments in recent years have shown that there is a number of hquids whose 
properties are highly anomalous [l|. For example, data for liquid metals, metalloids, non- 
metals, oxides and alloys, including Ga, Bi Te, S, Be, Mg, Ca, Sr, Ba, Si02, P, Se, Ce, Cs, 



observed in experiments on P, 
Data consistent with a L 



as water 



6|, silica 



a, 



c 



Rb, Co, Ge, GeisTegs 1|, reveal the presence of a temperature of maximum density (TMD) 
below which the density decreases under isobaric cooling. The most famous example of 
anomalous liquid is water, whose TMD at 1 atm is at approximately 4°C. Another anomaly 
that is attracting considerable attention is the possible existence of a liquid-liquid phase 
transition (LLPT) for single- component systems. Based on simulations for water [2|], it has 
been proposed that the two coexisting liquids, the high density liquid (HDL) and the low 
density liquid (LDL), would differ in density and local structure. Liquid polymorphism is 

31], P(OC6H5)3 |j] and yttrium oxide-aluminum melts |5|. 
!jPT arise from experiments for single-component systems, such 
8|, Se |9|, Co 10|, among others ll|. Here we ask the question 
how nano confinement of an anomalous liquid may change its behavior. The question is 
relevant for the vast arena of nanotechnological applications and for the lively debate about 
the properties of confined water in nanoscopic and biological systems 12l-ll5|. 

To approach the problem we focus on the effects of the structure of a fixed confining 
matrix of nanoscopic particles (NPs) on the phase diagram of a liquid with anomalous 
behavior and with a LLPT ending in a liquid-liquid critical point (LLCP). As model liquid 
we consider (i) a ramp potential IGj, and for comparison (ii) a shoulder potential 17l |. 

Potential (i) has a hard-core at distance r = a, and a linear ramp for a < r < 6 decreasing 
from interaction energy Ur > to f/o < 0, plus a linear ramp for b < r < c increasing from 
Ur to 0. We adopt b/a = 1.72, c/a = 3 and Ur = 3.56?7o. For this set of parameters the 
LLCP lies well above the melting line [l8|. The interaction with the NPs is modeled by the 
hard sphere repulsion. 

Potential (ii) has a repulsive shoulder and an attractive well with energy minimum Ur,, 

n n 

[171 ] with parameters chosen to fit a potential proposed for water \i^. The interaction with 
NPs is given by a l/r^^^ power law that crosses the zero of the interaction energy at the 
unit distance a. Results presented in the following, if not otherwise indicated, are for the 
potential (i), with consistent results for potential (ii). 



Both potentials (i) and (ii) display water-like anomalies and a LLPT 16|, [iTj. For both 
we perform simulations at constant number A^ of particles, constant volume V and constant 



temperature T, with periodic boundary conditions. For (i) we employ a discrete molecular 
dynamics algorithm by discretizing the linear ramp potential into steps, with Af/ = f/o/8 = 



0.125 20|. For (ii) we integrate the equations of motion using velocity Verlet integrator and 
we adopt the Allen thermostat 17|. 

We confine the liquid particles in a fixed matrix of A'np = 64 nanoparticles of diameter 
-Dnp = 3a. The liquid particles interact with NPs via hard core repulsion at distance 
To = [a + Dnp)/2. We consider three different structures for the matrix: a perfect cubic 
lattice (CUBE); a cubic lattice with Gaussian distortions (DIST) with a standard deviation 
of 1/4 of the separation between centers of NPs, which still preserves an approximately 
periodic and ordered structure of the confinement (FigH^); and a completely random (RND) 
configuration of NPs obtained by simulating a gas of hard spheres (Fig|T]D). The volume 
fraction of NPs is xnp = Vnp/^ = 24.5% in a cubic simulation box with fixed volume 
V/a^ = (20.6)^, and V^p/a^ = N^p47r{rQ/aY/3 volume of NP. We repeat the analysis for 
a random configuration of NPs with diameter 5a, keeping constant xnp and the size of the 
simulation box. We control p = N/{V — Vnp) of hquid particles by adjusting A^. We account 
for the excluded volume by rescaling the pressure P by V/iV — Vnp). 

For potential (i), the bulk hquid displays a LLCP at ^sT^^'^^/f/o = 0.375, P^^^^^aV f/o = 



0.243, and p^^^a^ = 0.37 [16|. Figures [T]^c) and [T]^d) show simulated isochores for DIST 
and RND confinement, respectively, with the HDL-LDL spinodal lines calculated from 
{dP/dp)T = and {d'^P/dp^)T 7^ 0, and the LLCP obtained from the merging of the 
spinodal lines where {dP/dp)T = {d'^P/dp^)T = 0. For each confinement we find that the 
LLCP shifts to lower T, higher p, and higher P relative to the bulk liquid (FigJJK). As the 
disorder in the confining matrix increases, T shift becomes more pronounced while p and 
P shifts are reduced. We find the same qualitative trend in the LLCP shifts for the liquid 
with potential (ii), and a progressive approach of the LLCP to the bulk case when the NPs 
concentration decreases (FigJJb), consistent with previous results for NP-liquid mixtures 

Although the periodic DIST confinement preserves the LDL-HDL coexistence region ob- 



served in bulk liquid 16|], which is consistent with a strong first-order LLPT, the RND 
confinement shrinks the coexistence region (Figs. [Th and[T]i), hinting at a weakening in the 
LLPT. This weakening is qualitatively consistent with recent numerical results for super- 



cooled water in random hydrophobic pore-like confinement 15 1. 



The region of density anomaly is bounded by the hnes of TMD and temperature of 
minimum density (TminD) located by the extrema of the isochores. In the bulk system the 
TminD line for high densities is hindered by the glass transition line and cannot be observed 
in the equilibrium liquid. Here we observe that the periodic structure of the confinement can 
dramatically affect density anomaly manifestations. Comparatively to the bulk, confinement 
decreases TMD and increases TminD, shrinking the T range of the density anomaly. While 
for the DIST confinement the density anomaly is still well defined, it appears much less 
pronounced in RND case. For a RND matrix of large confining NPs with diameter 5a, TMD 
and TminD are completely absent (not shown). 

To understand the origin of the different effects for the different confinements, we first 
study the density of the liquid in the vicinity of NPs. We find that a layer of liquid ad- 
sorbs onto the NPs, as emphasized by the NP-liquid particle radial distribution functions 
5'NP-iiq(^)(Fig|3K). We understand the increase of density near the NP surface as a conse- 
quence of the maximization of entropy. By sticking to the fixed NPs, the adsorbed liquid 
particles allow more free space to the the rest of the liquid, maximizing the entropy of the 
system (depletion effect). This result evokes a similar effect found for water forming a high- 
density layer of ~ 5 A at confining surfaces, regardless of the hydrophobic or hydrophilic 
interaction with the surface 



21 



23|. 



We find that, by increasing randomness in the confinement, the probability of overlap 
of NP exclusion volumes increases and the depletion effect decreases. As a consequence, 
the density of liquid near the NPs decreases (FigJSb). In addition, we analyze the density 
fluctuations and the associated measurable response function, the isothermal compressibility 
Krp (Figj3^), of the liquid in the vicinity of the NPs. We find that the local Kt is extremely 
small right at the interface, consistent with a tight packing of liquid particles around the 
NPs. Near the first minimum of 5'NP-iiq('") the local Kt is, instead, twice as high as the 



bulk, consistent with [2]J] . High vicinal density of liquid rationalizes the density increase of 
the LLCP (Figj2]), because, when part of the liquid is adsorbed onto the NPs, an average 
liquid density larger than bulk is necessary to build up the critical fluctuations. The shift 
is more pronounced for CUBE and DIST confinement, with respect to RND, because the 
more ordered the confinement the larger the NP surface available for the depletion effect. 

To understand the distribution of local density inside the simulation box, we performed 
Delaunay tessellation for a set of NPs in both RND and DIST confinements, partitioning 



the total volume into space-filling irregular tetrahedra with four nearest neighbor NPs as 
vertices (FigJ3^,b insets). We find that the RND confinement, with respect to the DIST case, 
induces a larger heterogeneity in the distribution of local densities p^ in each tetrahedra at 
any fixed value of the global density p (FigJ3^,b). The significant overlap of p„ distributions 
for the HDL and the LDL phases in RND confinement (FigJl^) suggests the presence of 
disordered mesoporous microdomains of various densities (FigiUi insets) that contribute 
to the weakening of the LLPT. This is not observed in DIST case (FigJDD), where the 
distributions of p„ are closely centered around the average global p with no overlap between 
HDL and LDL phases, preserving a strong LLPT. 

The fact that for RND confinement the local density is more heterogeneous and the 
liquid-liquid coexistence occurs at the local scale within Delaunay tetrahedra, allows us to 
rationalize the T-decrease of the LLCP with respect to the more ordered confinements and 
the bulk. Local heterogeneity inhibits the global LLPT, driving the LLCP temperature to 
a lower value. Due to the positive slope of the LLPT, the decrease of LLCP temperature 
implies a decrease of the LLCP pressure (Fig. [2]). The presence of density heterogeneity, 
along with the reduced depletion effect in the RND confinement matrix, gives us the key to 
understand the differences between ordered and disordered confinement. 

In conclusion, we predict that anomalous liquids, with a LLPT, preserve their bulk phase 
diagram and density anomalies when are confined in a matrix of periodically ordered NPs. 
Moreover, in the case of small distortion of the ordered confinement, the homogeneous crystal 
nucleation is inhibited at temperatures lower than in bulk, allowing the direct observation 
of the TminD locus. A strong depletion effect induces a large increase of density in the 
vicinity of the NPs. The effect is smaller when the confinement has a random structure. 
Randomness induces the occurrence of disordered density domains weakening the LLPT, 
narrowing the LLPT coexistence region, and washing out the density anomalies. 

It is interesting to mention that these results could qualitatively explain recent exper- 
iments for confined water, the prototypical anomalous liquid. While the TminD locus 
has been observed in supercooled water under hydrophilic confinement by MCM-41 silica 



nanoporous matrix [2J], its absence has been reported in hydrophobic mesoporous material 



CMK 



25|. MCM-41 forms a very regular matrix 2J], while CMK consists of grains, each 



with a disordered, 3D interconnected bicontinuous, pore structure 25|. Therefore, as shown 



in our model for a generic anomalous liquid, the disparity of results for different kinds of 



confinements may arise from the different amount of disorder in the confining structures. 
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FIG. 1: Effect of confinenient. Snapshots of the anomalous hquid (green) confined in a fixed 
matrix of NPs (yellow) in a DIST (a) and RND (b) configuration. Polynomial fits of simulated 
isochores of densities 0.89 < p/p^^^ < 1.59 (bottom to top at T/Tc = 1) for DIST (c) and RND (d). 
Randomness reduces the temperature and pressure of the LLCP (circles), the separation between 
the HDL (lower filled triangles) and LDL spinodals (upper open triangles) and the separation 
between the TMD (diamonds) and the temperature of minima of densities (TminD, squares). 
Samples of error bars on P are given in panel (c). Lines connecting symbols are guides for the 
eyes. Insets: cuts through the simulation box at density above (I), between (//) and below (///) 
the coexistence densities. The black circles (with a maximal radius of 3.5a) mark regions where 
liquid particles stick to the NPs, as discussed in the text. Outside these regions the liquid density 
is represented by a color code going from blue (high liquid density) to red (low liquid density). 
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FIG. 2: The confinement affects the parameters of the LLCP. Color-coded circles represent the 
LLCP parameters in the P — T — p phase space (a) for the liquid confined in the fixed matrix 
of NPs with CUBE, DIST and RND configuration. Increasing disorder in the confinement, from 
CUBE to DIST to RND, shifts the LLCP down in p, T, and P. (b) Upon decreasing concentration 
xnp (label near the symbols) for the CUBE confinement of the liquid with potential (ii), the LLCP 
approaches the bulk case. We find the same behavior for potential (i). 
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FIG. 3: (a) The liquid adsorbs onto the NPs. The NP-liquid particle radial distribution functions 
ffNP-iiq(^) at T/T^"-'^^ = 1.12 for CUBE (leftmost), DIST (center) and RND (rightmost) confine- 
ments for density p/ f^"^^ = 1.59 (thick lines) display large maxima (in parenthesis) at the closest 
NP-liquid particle distance r = r^ = 2a. Normalized local compressibility KxksT (broken lines) 
computed from the fluctuations of the number of particles in volume 5.5o'^ shows large peaks near 
the minimum of 5'NP-iiq(^)- The results for different confinements are shifted horizontally for clar- 
ity, (b) Number of liquid particles adsorbed onto a NP, at a distance tq < r < 2.5a, where the first 
minimum of (/NP-liq (?") occurs, for DIST (blue squares) and RND (red circles) for a range of densi- 
ties at fixed rescaled temperature r = (T^^D _ j^RNDyj^RND ^ (^^DIST _ j^piSTyj^piST ^ _qi^ 

Arrows mark the location of the critical density p^^^ (red) for RND and pj?^^""" (blue) for DIST. 
Lines are guides for the eyes. (Inset) The rescaled LDL-HDL coexistence region, with the chosen 
value for r, for DIST (gray) and RND (white) confinements, ending in the LLCP (circles). 
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FIG. 4: RND confinement induces heterogeneity in local densities, (a) The distribution of local 
density py/ f}^^^^ of the liquid included in Delaunay tetrahedra of the RND confinement for the 
LDL global density p/ f^^^^ = 0.91 (white histogram) partially overlaps with that for the HDL 
pjphnXk _ -^r^Q (black histogram), (b) The overlap is not present for the DIST case. In both 
panels an arrow marks the critical density, p^ for RND and pj? ^'^^ DIST. Insets: Delaunay 
tetrahedra for RND (top) and DIST (bottom) case, colored according to the value of p^, from high 
(dark gray) to low (light gray). 
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